#include <stdexcept>                              // std::runtime_error()

// STL
//#include <fstream>                              // std::ofstream
#include <functional>                             // std::function<>, ::bind(), ::placeholders::_x, ::cref()
#include <tuple>                                  // std::tie()
#include <vector>                                 // std::vector<>

// jScience
//#include "jmath_basic.h"                        // sign()
#include "jrandnum.hpp"                           // rand_num_uniform_Mersenne_twister()
//#include "jutility.hpp"                         // get_n_rand_unique_elements()

// stats++
#include "statsxx/dataset.hpp"                    // DataSet
#include "statsxx/machine_learning/NeuralNet.hpp" // NEURAL_NET, .()
#include "statsxx/optimization/EA.hpp"            // EAParam, EA()
//#include "statsxx/machine_learning/NeuralNet/utility.hpp" // fit_validation_error()


std::vector<double> new_individual(
                                   // -----
                                   NEURAL_NET mlp
                                   );

double new_gene(
                const int                  i,
                const std::vector<double> &w,
                // -----
                NEURAL_NET                 mlp
                );

double mutate_gene(
                   const int                  i,
                   const std::vector<double> &w,
                   // -----
                   NEURAL_NET                 mlp
                   );

double fitness_function(
                        const std::vector<double> &w,
                        // -----
                        NEURAL_NET                 mlp,
                        // -----
                        const DataSet             &ds_tr
                        );

bool convergence_function(
                          const std::vector<double> &w,
                          // -----
                          NEURAL_NET                 mlp,
                          // -----
                          const DataSet             &ds_val
                          );


//
// DESC:
//
// INPUT:
//
//          NEURAL_NET mlp    :: "template" neural network
//                               NOTE: This is used as a template for individuals, new genes and mutations, and ultimately to return.
//
//          EAParam    param  :: training parameters
//
//          DataSet    ds_tr  :: training dataset
//          DataSet    ds_val :: validation dataset
//
// OUTPUT:
//
//          NEURAL_NET mlp    :: (optimized) neural network
//
// -----
//
// NOTE: There is no "early stopping" for EA training.
//
// -----
//
// NOTE: This subroutine exists outside of NeuralNet, since it doesn't operate on a single object.
//
// -----
//
// TODO: NOTE: Perhaps mlp should be used to "initialize" new individuals, etc.
//
inline NEURAL_NET NeuralNet_EA_train(
                                     NEURAL_NET     mlp,
                                     // -----
                                     EAParam       &param,
                                     // -----
                                     const DataSet &ds_tr,
                                     const DataSet &ds_val
                                     )
{
//    NEURAL_NET mlp;

    //=========================================================
    // WRAP FUNCTIONS
    //=========================================================

    std::function<std::vector<double>()> new_individual_wrap = std::bind(
                                                                         new_individual,
                                                                         // -----
                                                                         std::cref(mlp)
                                                                         );

    std::function<double(const int, const std::vector<double> &)> new_gene_wrap = std::bind(
                                                                                            new_gene,
                                                                                            std::placeholders::_1,
                                                                                            std::placeholders::_2,
                                                                                            // -----
                                                                                            std::cref(mlp)
                                                                                            );

    std::function<double(const int, const std::vector<double> &)> mutate_gene_wrap = std::bind(
                                                                                               mutate_gene,
                                                                                               std::placeholders::_1,
                                                                                               std::placeholders::_2,
                                                                                               // -----
                                                                                               std::cref(mlp)
                                                                                               );

    std::function<double(const std::vector<double> &)> fitness_function_wrap = std::bind(
                                                                                         fitness_function,
                                                                                         std::placeholders::_1,
                                                                                         // -----
                                                                                         std::cref(mlp),
                                                                                         // -----
                                                                                         std::cref(ds_tr)
                                                                                         );

    std::function<bool(const std::vector<double> &)> convergence_function_wrap = std::bind(
                                                                                           convergence_function,
                                                                                           std::placeholders::_1,
                                                                                           // -----
                                                                                           std::cref(mlp),
                                                                                           // -----
                                                                                           std::cref(ds_val)
                                                                                           );

    //=========================================================
    // EA
    //=========================================================

    throw std::runtime_error("error in NeuralNet_EA_train(): EA (temporarily) unavailable");

    std::vector<double> w_opt;
/*
    std::vector<double> w_opt = EA(
                                   new_individual_wrap,
                                   new_gene_wrap,
                                   mutate_gene_wrap,
                                   fitness_function_wrap,
                                   convergence_function_wrap,
                                   // -----
                                   param
                                   );
*/
    //=========================================================
    // WRAP FUNCTIONS
    //=========================================================

    mlp.assign_weights(
                       w_opt
                       );

    return mlp;
}


//
// DESC: Create a new individual.
//
// NOTE: This uses the optimum link weights from creating a neural network.
//
inline std::vector<double> new_individual(
                                          // -----
                                          NEURAL_NET mlp
                                          )
{
    mlp.optimize_link_weights();

    return mlp.get_weights();
}


//
// DESC: Create a new gene.
//
// NOTE: This uses the optimum link weights from creating a neural network.
//
// TODO: NOTE: ... This is fairly inefficient, as ALL weights must be reinitialized.
//
// TODO: NOTE: ... However, it keeps all weight initialization consistent, ...
// TODO: NOTE: ... and all functionality contained in NeuralNet.
//
inline double new_gene(
                       const int                  i,
                       const std::vector<double> &w,
                       // -----
                       NEURAL_NET                 mlp
                      )
{
    mlp.optimize_link_weights();

    return mlp.get_weights()[i];
}


//
// DESC: Mutate a gene.
//
// NOTE: This uses the optimum link weights from creating a neural network.
// NOTE: ... (For reasons discussed in the code below).
//
// TODO: NOTE: ... This is fairly inefficient, for reasons discussed above for new_gene().
// TODO: NOTE: ... But the same NOTEs apply.
//
inline double mutate_gene(
                          const int                  i,
                          const std::vector<double> &w,
                          // -----
                          NEURAL_NET                 mlp
                          )
{
    // NOTE: Keeping in mind that initialization of weights targets a mean of 0 and a uniform distribution dependent on how many weights are attached to the associated neuron, ...
    //
    // NOTE: ... a reasonable approach to mutation is (seems to be) to take the current weight and shift it by the same (initialization) amount.

    mlp.optimize_link_weights();

    return ( w[i] + mlp.get_weights()[i] );
}


//
// DESC: Fitness function.
//
inline double fitness_function(
                               const std::vector<double> &w,
                               // -----
                               NEURAL_NET                 mlp,
                               // -----
                               const DataSet             &ds_tr
                               )
{
    mlp.assign_weights(
                       w
                       );

    double err_mean;
    double err_stddev;
    std::vector<std::vector<double>> out;
    std::tie(
             err_mean,
             err_stddev,
             out
             ) = mlp.parse_TS(
                              ds_tr
                              );

    // NOTE: The following fitness maximum from the error minimum has the nice property that it is bounded in [0,1].
    //
    // TODO: NOTE: ... There is possible room for improvement here.

    return (1./(1. + err_mean));
}


//
// DESC: Convergence function.
//
inline bool convergence_function(
                                 const std::vector<double> &w,
                                 // -----
                                 NEURAL_NET                 mlp,
                                 // -----
                                 const DataSet             &ds_val
                                 )
{
    // TODO: Need to figure out.

    return false;
}
